MASSCLEAN - MASSive CLuster Evolution and ANalysis Package - 

Description and Tests 



Bogdan Popescu and M. M. Hanson 
University of Cincinnati, Physics Department, 400 G/P, PO Box 210011, Cincinnati, OH 45221-0011 
popes cbOemail .uc . edu, margaret . hansonOuc . edu 



ABSTRACT 

We present MASSCLEAN, a new, sophisticated and robust stellar cluster image and photome- 
try simulation package. This visualization tool is able to create color-magnitude diagrams and 
standard FITS images in any of the traditional optical and near-infrared bands based on cluster 
characteristics input by the user, including but not limited to distance, age, mass, radius and 
extinction. At the limit of very distant, unresolved clusters, we have checked the integrated colors 
created in MASSCLEAN against those from other simple stellar population models with consistent 
results. We have also tested models which provide a reasonable estimate of the field star con- 
tamination in images and color-magnitude diagrams. We demonstrate the package by simulating 
images and color-magnitude diagrams of well known massive Milky Way clusters and compare 
their appearance to real data. Because the algorithm populates the cluster with a discrete num- 
ber of tenable stars, it can be used as part of a Monte Carlo Method to derive the probabilistic 
range of characteristics (integrated colors, for example) consistent with a given cluster mass and 
age. Our simulation package is available for download and will run on any standard desktop 
running UNIX/Linux. Full documentation on installation and its use is also available. Finally, a 
web-based version of MASSCLEAN which can be immediately used and is sufficiently adaptable for 
most applications is available through a web interface. 
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1. Introduction 

Stellar clusters provide among the most critical 
observational benchmarks for testing the physics 
of stellar and galactic evolution and galactic struc- 
ture. Stellar clusters are often exploited by as- 
suming their constituent stars all formed from the 
same interstellar material (similar initi al chem 



i stry) and at the same time (same age) . iTinslev 
(|1978l ) pioneered the application of using stellar 
evolutionary codes to model the observed charac- 
teristics of unresolved stellar clusters and galax- 
ies, referred to as evolution ary population synthe- 
sis (jTinslev fc Gunnlll976T ). The present day ver- 
sion of such methods is seen in the application of 
simple stellar population (hereafter SSP) models 



to unresolved systems. 

Owing to their great utility, a vast array of 
SSP models have been dev eloped by researchers 

in the past two decades: Leitherer et al.l 19991: 

Hurley et al.l2000HSchulz et al.120021: iBruzual fc Charlotl 



20031: IVazquez fc Leithereil l2005t iMaraston 1 120051 



to name just a few of the most widely used. 
As these models are applied to a multitude of 
galactic and extragalactic applications, they each 
have slightly different techniques in how they in- 
clude the input physics (beyond differences in 
which stellar evolutionary codes they choose). 
For example, the most common method derives 
the emergent properties by integrating along a 
single stellar isochrone - isochrone synthesis - 
(jCharlot fc Bruzuall 1 1 9 9 lh . while another tech- 
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nique follows the cluster evolu tion based on fuel 
consumption ( Maraston 2005h . This can lead to 



perceivable differences in predicted outputs for 
the same input cluster characteristic (age, metal- 
licity, stellar mass function) from one SSP code 
to the next, even when the same stellar evolu- 
tionary codes are used (see the tests p resented in 
Beaslev et aT1l2002l IPessev et al.ll2008l ). 

Among the variety of SSP codes already avail- 
able, we are unaware of any that take advan- 
tage of the situation when individual stars are 
fully or at least partially resolved. For exam- 
ple, the GALEV evolutionary sy nthesis mod- 
els o f Fritze and collaborators (ISchulz et al 



2002t [Anders fc Fritze - v. Alvenslebenl 12003: 



Kotulla et al.1 120091 ) can produce very useful syn- 
thetic CMDs for discrete stellar clusters (see Fig. 



A3 in iKotulla et alj 120091 ) though they have not 



gone so far as to attempt to make simulated im- 
ages of stellar clusters. The need for a stellar clus- 
ter image simulation, based on the tenets of the 
traditional SSP code, drove us to develop a new 
analysis tool. Our motivation to develop an im- 
age simulation code was to apply it to the search 
and analysis of deeply embedded massive open 
clusters lying in the inner Milky Way. However, 
the code can be applied to Local Group galaxies 
where massive clusters can be partially resolved 
with some telescopes. 

Our visualization and analysis tool, called 
MASSCLEAN, provides image simulations and thus 
can be used to answer entirely different, important 
questions in stellar astronomy and galactic struc- 
ture than previous SSP codes, both within the 
Galaxy and in nearby external galaxies. In this 
paper we present the details and current coding 
of the first release of our new, SSP-likc imaging 
and photometric simulation code, MASSCLEAN. In 
SJ21 we provide a quantitative description of the 
computational algorithms used in the routine. In 
[J3]we provide test runs of the simulations compar- 
ing the integrated colors and magnitudes against 
those of currently used SSP codes. We also provide 
example images with an eye towards accurately 
modeling well known massive young clusters. In 
[]H we discuss both the promise and limitations of 
MASSCLEAN for a variety of astrophysical applica- 
tions. 



2. General description of MASSCLEAN 

Our simulation package, MASSCLEAN^, uses a 
nominal number of input parameters: mass, initial 
mass function (hereafter, IMF), metallicity, ex- 
tinction, distance, spatial distribution parameters 
and stellar field density. The predicted character- 
istics are computed for a range of cluster ages from 
10 6 to a few 10 9 years. The user can also choose 
to include a parameter that allows for mass seg- 
regation as the cluster ages. Many other features 
are described in the next sections. 

The simulation code is built using numer- 
ous well established theoretical and empirical 
models for stars and stellar clusters, begin- 
ning with the Kr oupa-Salpeter IMF for stella r 
mass distribution ( Kroupa 2002 : Salpeter 1955t ). 
The user has the option to choose between two 
stellar evolution models: the Geneva Models 



( Leieune fc Schacrcr 2001) or the Padova Models 
( Marigo et all2008l). The extinct i on m odel based 
on ICardelli. Clayton and Mathisl I1989L t he King 



Model for spatial distribution (King 1962) an d the 
SKY Model for the stellar field dWainscoat et al ' 



Il992t ICohenlll994t ICohenl [T995: Bahcall fc Soneira 



119841 ) are further used. This is reviewed in more 
detail below. 

It was our intent to make MASSCLEAN user 
friendly and versatile. The package is designed 
to allow the user significant latitude and flexi- 
bility in how they will apply it. MASSCLEAN is 
written as a series of independently run sub- 
routines, performing individual calculations at 
various stages. This allows users to substitute 
their own inputs at any stage or to skip or per- 
form their own calculations in lieu of those pro- 
vided within the MASSCLEAN package. Information 
is passed from one routine to the next through 
read/write of ASCII files. While these files be- 
come quite large (gigabyte is not uncommon for 
the entire output of hundred of files), this allows 
users of the MASSCLEAN routines to easily check, 
edit and or substitute the output at any stage in 
the calculations. Although it can take some time 
to generate images, the code for the photomet- 
ric simulation is very fast. The package is freely 
available under GNU General Public License at 
http : //www. physics .uc . edu/~bogdan/. Down- 
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loadable documentation on installing and running 
the code is also available. Finally, MASSCLEAN. web_, 
a web-based interface, is immediately available^ 
and can be used for many basic applications. 

2.1. The Mass Distribution 

The number of stars formed in the (M, M+dM) 
range is: 



dN = £(M)dM 



(1) 



The multi-part po wer law £(M) derived from 
Kroupa-Salpeter IMF (jKroupa 2002HSalpeter 19551) 



is : 



£(M) = k < 



(^y ai > m <M< mi 

(^y 02 > m!<M<m 2 

Imy^i^y 03 > -2<M<m 3 

(2) 

with mass expressed in Mq units and : 

a x = +0.3 ± 0.7 , 0.01 < M/Mq < 0.08 
a 2 = +1.3 ± 0.5 , 0.08 < M/M© < 0.50 
"3 =il?lo 3 . 0.50 < M/M < m 3 



and m3 = oo or 



(3) 

(for an I MF with 



upper mass cutoff (lOev fc Clarke 2005h ). Using 
{{M)/k = £i(M) (with i = 1,2,3 respectively), 
the total mass of the cluster can be written : 



M, 



total 



M(N)dN 



(4) 



M total = I M—AM = / £(M)MdM 



,<;<;, CIA/ . / 



(5) 



M totQ , =J] fe/ 6(M)MdM (6) 

i=l \ •'mi- 1 / 



The normalization constant 



k = 



Eli (.C ; C':.W :.U,I.U 



(7) 



2 http: / /www .physics .uc . edu/~bogdan/masscleaii/ 



From the equations (JTJ) and ([7]) 



TVi (M, M + dM) 



total J M 



Ci(M)dM 



TLi (€lM M ) MAM 



(8) 

Our package uses equation (JSJ) to compute the 
mass distribution of stars based on the total mass 
of the cluster M to tai and the IMF, described by 
the three £j(M). The mass range to be included 
in the distribution is chosen by the user. The val- 
ues on and rrii are also input parameters. In this 
way, a Kroupa IMF can become a Salpeter IMF 
when ot\ = a 2 = Q3 = 2.35. The versatility of the 
package allows us to use a first-order, second-order 
or third-order power law model (Kroupa-Salpeter 
type) for the IMF. The use of a truncated IMF 
is optional. All the mass bins (M, M + dM) are 
computed such that the value of Ni(M, M + dM) 
is an integer (up to some tolerance, also chosen by 
the user). The program allows fluctuations in the 
computed mass distribution, such that the distri- 
bution is different for every run. The user can turn 
this feature on or off. 

The characteristics of the stellar mass distribu- 
tion, whether truncated or tied to the cluster mass, 
etc., are an important component of a stellar clus- 
ter for a variety of investigations. This is why 
the MASSCLEAN package offers a wide array of op- 
tions to build the stellar mass distribution. Sam- 
ple computed mass distributions are presented in 
Figure Q] 

2.2. Evolutionary Models 

A variety of evolutionary tracks have been used 
in simulation codes to define the relationship be- 
tween age and metallicity and int egrated broad- 



1996, 



band colors o r spec tr al features dCharlod 
Brocato et ail 120001 iBruzual k Chariot 
Leitherer et al.lll999l ). Selecting which models to 



2003, 



use is influenced by the goals of that simulation. 
Presently, the MASSCLEAN package supports two 
different evolutionary models and the user can add 
their own. In this first version of the MASSCLEAN 
package, we have chosen to inclu de the Geneva 
Models (jLejeune fc Schaererl 120011 1 and also the 



recen tly released Padova Models (jMarigo et al 



2008) The Geneva models provide excellent treat- 
ment of the evolutionary properties of high mass 
stars and ample time sampling to closely follow 
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M (solar masses) 

(a) ai = 0.3, Q2 = 1-3, as = 2.3 



M (solar masses) 

(b) qi = 0.3, a 2 = 1.3, a 3 = 2.1 



Fig. 1. — Mass distribution for M to tai = 5 x 1O 4 M 
using Kroupa IMF. Because all the bins above 7Mq 
contain only one star, only the [0.8, 7.0] interval is dis- 
played. The solid line shows the mathematical form 
of the distribution for constant width bins. 



the evolution of a young cluster (age > 10 7 years). 
This is important in the study of young clusters 
where evolution occurs quickly and the light from 
a few massive stars can dominate the light of the 
cluster. The Padova models are considered supe- 
rior in their treatment of older stars, as it more 
carefully considers the evolution of AGB popu- 
lations and several peculiarities present at the 
onset of the thermally-pulsating asymptotic giant 
branch phase (TP- AGB). The package can easily 
be switched to use a different set of isochrones 
and tracks as the user sees fit, such as evolution- 
ary models taylored to old stars, or low metallic- 
ity, etc. One can also expand the current set of 
photometric bands to include the ultraviolet and 
mid-infrared. 

2.3. CCM Extinction Model 

We wish to use our simulations in the study 
of Milky Way clusters as well as more distant 



extragalactic clusters. In both applications, we 
must pay close attention to extinction effects. 
For the extinction, the user can enter the value 
in two ways: the exact extinction value can 
be entered manually (in the configuration file) 
for every band or the extinction in each band 
can be computed using an inputted Rv and 
Ay and applying the CCM Extin ction Model 
(jCardelli. Clayton and Mathid . Il989h . A sample 
computed extinction curve is presented in Figure 

m 
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X (nm) 

Fig. 2.— CCM Extinction Model, A v = 11.6, Rv 
3.1 



2.4. Spatial Distribution 

The MASSCLEAN package is unique to previous 
SSP codes in that it produces simulated images of 
the stellar clusters. Thus, the code is concerned 
with selecting appropriate spatial parameters de- 
fined by the angular size (linear scale and distance 
to the cluster) and for the first time, intrinsic stel- 
lar density. To accomplish this we have intr oduced 
the King Model Distribution (jKing 1962T ) . given 
by: 





f(r) = k 



An anisotropic spatial distribution can even be 
generated such that it still obeys the radial King 
profile. An ellipsoid in a prolate or oblate projec- 
tion can be generated according to the parameters 
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entered by the user in the configuration file. The 
rotation angle is also selected by the user. 

A simple linear mass segregation feature allows 
more massive stars to fall toward the cluster's cen- 
ter as the cluster ages (this feature can be turned 
on or off). 

The file containing the spatial distribution can 
be replaced by another file provided by the user, 
for example, based on coordinates from a real im- 
age. 

Sample computed spatial distributions are pre- 
sented in Figure [3] 



(a) Isotropic Distribution (b) Anisotropic Distribution 



Fig. 3.— King Distribution, r c = 0.50 pc, r t = 0.72 



The second option is to use a real stellar field. 
Since the file format is the same as the one used by 
SkyMaker (lBertinll200il: iBertin fc Fouquell2007t ) 
and SExtractor ( Bertin fc ArnoutHl996f) . from 
the computational point of view there is no differ- 
ence between a real and a simulated stellar field. 
MASSCLEAN can use a real stellar field provided by 
SExtractor. The user can also choose not to 
include any stellar field at all. 



22 21 20 19 18 17 16 15 14 13 12 11 10 9 8 
mag 

Fig. 4. — Cumulative number of stars vs. magnitude 



pc, Mtotai = 10 J M G 



2.5. Description of the Stellar Field 

The MASSCLEAN package has the option to in- 
clude a field star population. This may be a sim- 
ulated one, using the SKY Model for the stellar 



field i|Wainscoat et alJll992l ICohenl Il994t ICohen 
1995t Bahcall fc SoneirJ 1984) or a real one. 



In the first case, starting from the total num- 
ber of stars brighter than the selected magnitude 
limit (which is an input parameter), the distri- 
bution is generated using the slope of the cu- 
mulative numbers of stars (which is also an in- 
put parameter), such as shown in Figure |H Col- 
ors are computed based o n the Geneva Models 
dLeieune fc Schaererl l200ll) and the BaSeL~2.2 
grid (Leieune et al.l Il998l ). and Padova Models 
(iMarigo et al.ll2008l) and the ATLAS9-ODFNEW 
grid ( Castelli fc Kuruc^ 12003 ). respectivelly. Ex- 
tinction is also included and the user has the op- 
tion to provide the values for every band. The 
program used to generate Figure 0] can compute 
the necessary parameters for a simulated stellar 
field based on the properties derived from a real 



2.6. Creating the FITS Images and HR 
Diagrams 

Using all of the above models, MASSCLEAN com- 
putes actual mass, absolute and apparent mag- 
nitude (UBVRIJHK), color indexes, temperature, 
luminosity and position for all the stars and all the 
ages included in the evolutionary database (which 
ever isochrones one employs). The default age 
range is log(age/yr) = 6.0 — 9.5, but the option to 
run it only for a few selected values is also avail- 
able. These outputs are directly used to generate 
the color-magnitude diagrams and images. There 
are available relations of transformation from UB- 
VRIJHK to different bandpasses for the Geneva 
and Padova models. It is expected that the band- 
passes for HST, 2MASS and Spitzer will be di- 
rectly available soon (and will be included in a 
future upgrade). 

The FITS images in each of the broad bands 



are generated using SkyMaker ijBertinl 12001 



Bertin fc Fouque1l2007i ). MASSCLEAN writes all the 
necessary scripts to run SkyMaker that gener- 
ates the images. The configuration files for galac- 
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tic and extragalactic clusters are also provided, so 
no knowledge about SkyMaker is required. The 
description of the point spread function (PSF) is 
available in these configuration files and can be 
changed by the user. In addition to that, Sky- 
Maker can work with a different PSF provided by 
the user in a separate file. The users can make any 
plot from the package's output using the plotting 
program of their choice, but MASSCLEAN also writes 
a script which generates HR and color-magnitude 
diagrams. 

3. Tests of the MASSCLEAN Package 

Before demonstrating the unique utility of 
MASSCLEAN and directly comparing its output 
to real clusters, we will first provide a few log- 
ical tests. Our first test will be to derive from 
MASSCLEAN, the same kind of values which come 
from other stellar population simulation programs 
which are widely used in the field to ensure that 
MASSCLEAN gives consistent results. MASSCLEAN 
does not generate sp ectral features, such as pro- 
vided by Starburst99 (jLeitherer et all 19991 ). How- 
ever, it can be made to generate integrated mag- 
nitudes and colors. This is achieved by simply 
summing up the flux over all stars in the final 
simulated cluster. 

3.1. Integrated colors as a function of age 
and metallicity 

The most commonly-used discriminator in 
the study of extragalacti c super cluster stud 



ies is integrated col ors (jHoltzman et al.1 |1992, 



Fusi Becci et al.l 120051 ). This is because it is typ 
ically the only information ob tainable for very dis- 
tant, unresolved star clusters (jWhitmore k, Schweizer 
1995). Among the most widely-used SSP mod- 
els for i nterpreting integra t ed co lors are those 
given bvlBruzual fc Chariot) d2003h ( GALAXEV), 



Maraston 



. .20051 GALEV 2003 (ISchulz et all 
2002l and l Anders & Fritze - v Alvenslebenll2003l) . 
GAL EV 2009 dKotulla et all l2009h. and Padova 
2008 (|Marigo et al.ll2008l ). We shall provide a side 



by side comparison of MASSCLEAN against these 
three modern SSP codes. 

The default age range for MASSCLEAN is [6.0, 9.5] 
in the logarithmic scale. In Figures [5] though [9] we 
use a smaller age range in order to accomodate the 
data from all the available models. In Figure[5]we 



plot the variation of integrated colors with age for 
different SSP models. Models using solar metallic- 
ity are presented here. In Figure [5] (a) the results 
from MASSCLEAN using Padova 2008 is shown as the 
solid line and compared with GALEV 2009 Mod- 
els (Padova) in the dotted line, GALAXEV mod- 
els (Padova 1994) are shown in the dashed-double 
dotted line, Maraston (2005) is shown in the dou- 
ble dotted line and Padova 2008 (Marigo et al. 
2008) is shown in a dashed line. Specifically, the 
MASSCLEAN colors follow quite closely to the col- 
ors given by GALAXEV (Padova) and are essen- 
tially identical to Padova (2008), as it should be. 
In Figure [5] (b) the results from MASSCLEAN using 
the Geneva evolutionary models is shown as the 
solid line and compared with GALEV 2003 Mod- 
els (Geneva) in dotted line, GALAXEV (Geneva) 
are shown in the dashed line, and Maraston (2005) 
is again shown in the double dotted line. Figured] 
(a) and (b) demonstrates that our code, based on a 
finite stellar generation algorithm, gives the same 
integrated color results of other SSP codes that use 
the same evolutionary models but using a statisti- 
cally weighted mass distribution. The MASSCLEAN 
integrated colors were computed in the high mass 
limit (~ 10 6 Mq), which allows for isochrones to 
be fully populated for all masses. 

In Figure [6] is presented the variation of inte- 
grated colors with age computed by MASSCLEAN 
for a moderate cluster mass (6.5 x 10 4 M©). Fig- 
ure[5](a) is using Padova 2008, and (b) is using the 
Geneva models. This is a demonstration the other 
SSP codes are not able to make since total mass 
is not available as an input. The bars represent 
the range observed in the output integrated colors 
from 1000 randomly generated cluster simulations 
(Monte Carlo runs of the stellar mass function). 
What is shown here is an increase in observed 
deviation in color with longer wavelength colors. 
This is further explored in Figure [7] Here, we've 
concentrated on the mass effect seen for the V- 
K colors only, as this band is where the observed 
color range is known to be greatest. These data 
in this figure were constructed by again creating 
1000 randomly generated clusters, now over a va- 
riety of masses. The figure provides the results 
from three general mass ranges. As cluster mass 
decreases, the observed range in integrated colors 
increases. 

Finally, the effect of changing the metallicity in 
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(a) (b) 

Fig. 5. — Integrated colors for different SSP Models as a function of age for solar metallicity: (a) Padova 2008 
(MASSCLEAN) - solid line; GALEV (Kotulla et al. 2009) - dotted line; GALAXEV (Bruzual & Chariot 2003) with 
Padova 1994 - dashed-double dotted line; Maraston (2005) - double dotted line; and Padova 2008, Mariago et al. 
2008 is the dashed line, (b) Geneva (MASSCLEAN) - solid line; GALAXEV (Bruzual & Chariot 2003) with Geneva - 
dashed line; GALEV (Anders & Fritze 2003) - dotted line; Maraston (2005) - double dotted line. 



the stellar evolution models (in (a) for Padova, in 
(b) for Geneva) is also displayed in Figure [7] The 
variation of integrated colors for different IMFs as 
computed by MASSCLEAN are presented in Figured] 
Previous studies have noted the very weak sensi- 
tivity of the IMF on the integrated colors of stellar 
clusters over this age range (see, f or example, Fig- 
ure 4 of iBruzual fc Charlotlfiool . 

In Figure [9] we show a full scale diagram 
of the integrated broadband colors, U-B, B-V, 
and V-K derived with MASSCLEAN models using 
Padova and Geneva, and for two metallicities. 
Here we also show observed data, taken from 



Hunter et al. (|2003l) . for stellar clusters in the 
Magellanic Clouds in the upper two panels. In 
the lower panel, we have borrowed the observed 
photometry given in iMarigo et al. ( 20081 ). Fig. 9, 
to show the observed variation in V-K color for 
LMC clusters. The most obvious deviations are 
seen with the prediction of fairly red clusters in 
the age range from log(age) = 8 to 8.8 by the 
Pa dova evolutionary models (this was first noted 
by[M 

arigo et al.ll2008l ) and the prediction of fairly 
blue clusters with log(age) > 9.0 by the Geneva 
models. 

The color dispersion presented in Figures [6] and 
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along with the photometric data presented in 
Figure [H show the advantage of the MASSCLEAN 
approach with respect to real clusters. As de- 
scribed by iMarigo et al. (|2008h (§5.4 and Fi gurc 
9), the continuous and smooth SSP lines are not 
able to describe accurately the observational data. 
A more accurate description of the color dispersion 
as a function of the mass of the cluster will be pre- 
sented in a subsequent paper (Popescu & Hanson, 
in prep.). 

Finally, in the top panel of Figure \W\ we present 
the evolution of integrated colors (U — B)q vs 
(B — V)o for several different SSP models along 
side predictions coming from MASSCLEAN (solid 
line). All simulations were made using solar metal- 
licity. The MASSCLEAN colors compare well to 



real data as shown already with previous work 
(jGirardi et al.lll995| ). The lower panel shows pre- 



dicted color evolution using MASSCLEAN, but show- 
ing a variety of metallicities. As in Figure [H the 
solid dots in the lower two panels represent ob- 
served photometric dat a for stellar clusters in the 
Magellanic Clouds from lHunter et a l. (2003). 



3.2. Color-magnitude diagrams 

Our next test will be to simulate HR and color- 
magnitude diagrams of clusters using MASSCLEAN. 
As a first demonstration, we provide in Figure QTJ 
an HR Diagram for a stellar cluster with a total 
mass of 10 5 M and a log(agefyr) — 6.65. While 
the individual stars in the cluster are clearly seen 
at the high mass end, crowding prevents one from 
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(a) (b) 

Fig. 7. — The range of variation in V-K for three mass intervals: 50 — 75-K" (bottom panel), 75 — 125K (middle panel), 
125 - 175K (upper panel), (a) Padova 2008 (MASSCLEAN) (b) Geneva (MASSCLEAN). The height of the bar shows the 
observed range of output colors from 1000 Monte Carlo runs where input mass was allowed to vary over the range 
listed. 



seeing anything but a broad blur of stars below 
about 40 Mq. Naturally, this figure looks no dif- 
ferent from the Geneva isochrones used for that 
same age. In Figure ITTT a) the stars have been 
rebinned to take on a more natural range of val- 
ues; in Figure fTTTb) the stars follow a pure single 
isochrone. We also provide an example of a color- 
magnitude diagram, V versus B-V for this same 
simulated cluster, in Figure EEH At the low mass 
end, again a small variation has been introduced 
to give the main sequence some width in Figure 
[T27 a). The widening in Figure fTTT a) and Figure 
I12f a) is a user's choice, the maximum width is an 
input parameter (this option can also be turned 
off) and could be set up to correspond to the pho- 



tometric error. 

The cluster mass and age shown in these two 
figures were selected rather explicitly to rep- 
resent a known cluster. It matches the esti- 
mated mass and age of the Milky Way clus- 
ter, Westerlund 1. The para meters used for th e 
CMD simulation are based on We sterlund! 
Clark fc Neeueruelal (l2002h . IClark et all <|2005h 
and iFiger et alj ( |2006t h We adopted the follow- 
ing values: M to tai = 10 5 M©, solar metallicity, 
Ay = 11.6 mag, Ry = 3.1, r t = 0.72 pc and 
r c = 0.50 pc. 

There is some uncertainty in the distance to 
Westerlund 1. We adopted a value of d = 4 kpc 
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Fig. 9. — MASSCLEAN integrated colors for differing evolutionary models and metallicities displayed with observed 
photometric data taken in the Magellanic Clouds from lHunter et all l|2003h . the solid dots, and lMarigo et all (2008), 
the solar symbols. 
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Fig. 10. — Integrated colors (U — B)o vs (B — V)o- The circles along the solar metallicity path in the lower panels 
correspond to the Alog(age/yr) = 0.05. 
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Fig. 8. — Integrated colors for different IMFs as a 
function of age computed by MASSCLEAN (Geneva): 
Kroupa IMF a 3 = 2.6 - solid line; Kroupa IMF 
«3 = 2.3 - dashed line; Salpeter IMF - dotted line 
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Fig. 11 — MASSCLEAN simulated HR Diagram, for 
cluster properties log(age) = 6.65, M to tai = 1O 5 M0. 
In (a), a rebinning has been done to give the stars a 
more realistic range of values. In (b), one sees the 
stars following a theoretical track from the Geneva 
Database. Both figures cover the stellar masses from 
120M© to 1M Q . 
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Fig. 12.— The V vs (B - 

gram, Av = 11.6, log(age) = 6.65, M to tai = 10 D M©. 
The tracks correspond to the same values of mass as 
in Figure [TT1 



(distance modulus 13.01 mag). A Kroupa IMF 
has been used, with ai = 0.3, a.2 = 1.3 and 0:3 = 
2.3. Our simulation shows that the best agreement 
with the actual data corresponds to the isochrone 
log(age) = 6.65. This value als o agre es with the 



recent results of iBrandner et al.l (|2008r ). 



In Figure [13] and Figure [14] we present near- 
infrared color-magnitude diagrams for Westerlund 
1. MASSCLEAN isochrones are overlayed with real 
photometry for the cluster. We cannot present 
optical color-magnitude diagrams for this clus- 
ter because of the enormous incompleteness due 
to high extinction (Ay = 11). The photometric 
data in Figure [T3] and Fig ure [T4] was taken from 



2004) 



the 2MASS and NOMAD (jZacharias et al. 
catalogs. 

If MASSCLEAN is to be used to constrain the 
properties of observed clusters, a method for deter- 
mining the best simulation inputs must still be de- 
veloped. There is unlikely to be a straightforward 
way to determine the best fitting simulated image. 
Rather, we expect comparisons between simulated 
and real data will need to be done within the 
CMDs, comparing photometry. This is presently 
being worked on, but we expect to base such a 
goodness of fit on observed versus simulated cumu- 
lative distribution functions in various photomet- 
ric bands and using a Kolmogorov-Smirnov test 
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to select the closest fitting simulation (Popescu & 
Hanson, in prep.). 
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Fig. 13. — H vs. J — H color-magnitude diagram 
for Westerlund 1. Circles - 2MASS source catalog 4.8' 
radius; dots - 2MASS source catalog 2.0' radius; tri- 
angles - NOMAD catalog 4.8' radius. The solid lines 
correspond to the Geneva isochrones (log(age/yr) = 
6.65,6.75), and the dotted lines correspond to the 
Padova isochrones (log(agefyr) — 6.70,6.75). 
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Fig. 14. — K vs. H — K color-magnitude diagram 
for Westerlund 1. Circles - 2MASS source catalog 4.8' 
radius; dots - 2MASS source catalog 2.0' radius; tri- 
angles - NOMAD catalog 4.8' radius. The solid lines 
correspond to the Geneva isochrones (log(age/yr) = 
6.65,6.75), and the dotted lines correspond to the 
Padova isochrones (log(agefyr) — 6.70,6.75). 



3.3. Image simulations of Galactic clusters 

Our most visual demonstration will be the 
image simulations provided by MASSCLEAN. The 
MASSCLEAN package has been used to simulate sev- 
eral well known young Milky Way clusters: NGC 
3603, h and X Persei {NGC 869, NGC 884), as 
well as Westerlund 1. This corresponds to a mass 
range from 4.3 x 10 3 M Q to 1O 5 M . 

Using the input characteristics given in §3.2| 
we have created a simulation of Westerlund 
1 in the J-band and provided in Figure [T5] 
(a). This can be compared directly with the 
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2MAS^I( Skrutskie et alJ l2006h J -band image of 
Westerlund 1 shown as Figure [15] (b) . In look- 
ing at this comparison, one might also consider 
Figs. [TTTa) and [T27 a). Presently, the simulated 
image has a non-realistic 'flatness' in the stellar 
brightnesses due to several RSGs exhibiting near 
identical magnitude and color. Some amount of 
rebinning to randomize slightly the properties of 
the RSGs (something the cluster will have done 
more natually with the non-zero cluster age dis- 
tribution) may be desirable to include in future 
updates and releases of the code. 

Using the data from IStolte et all (|2004l l2006h 
we simulated NGC 3603 using: log(age) = 6.00, 
Mtotai = 10 4 M Q and solar metallicity, Ay = 4.5 
mag, Rv — 3.1, distance d — 6 kpc (distance 
modulus 13.9 mag). For the spatial distribution 
we used rt = 4.4' and r c — 0.4'. We note that 
a single-power Salpeter IMF could not lead to 
an agreement with the actual photometric data. 
The best fit corresponds to Kroupa IMF with 
a\ = 0.3, a2 = 1.3 and 0:3 = 2.4. The J-band 
i mage is shown in Fig ure [TBT a). and the 2MASS 
( Skrutskie et al . 200G) image is shown in Figure 
IToT b). It should be clear to the reader we do not 
include nebulosity in the image simulations. 

We have built a suite of different cluster mod- 
els by varying all of the input variables within 
the range of known measurements for these two 
clusters. These resulting simulations were tested 
against the available catalogs of images and photo- 
metric data for the two clusters as a means of con- 
straining the cluster properties. The phase space 
of properties is far too large to enable MASSCLEAN 
to serve as an efficient means of modeling clus- 
ters. However, if a reasonable guess can be made 
as to many of the cluster properties (age, distance, 
extinction) MASSCLEAN can work extremely effec- 
tively in a limited region of characteristics to de- 
termine the best model match to the observed 
data. A complete analysis for NGC 3603 and 
Westerlund 1 will be presented in a subsequent 
paper, along with the development of a proper 
goodness of fit (Popescu & Hanson, in prep.). 



3 This publication makes use of data products from the Two 
Micron All Sky Survey, which is a joint project of the Uni- 
versity of Massachusetts and the Infrared Processing and 
Analysis Center/California Institute of Technology, funded 
by the National Aeronautics and Space Administration and 
the National Science Foundation. 



(a) Simulated Image 



W m 




(b) 2MASS Image 

Fig. 15. — MASSCLEAN simulated images of Wester- 
lund 1. In the upper panel (a) we provide the sim- 
ulated image in the J Band. A simulated field star 
population has been included. The lower panel comes 
from 2MASS and is also in the J Band. Both fields 
shown are for a 4.8' x 4.8' image. 
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3.4. Simulations of multiple starburst 
events or subclustering 
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MASSCLEAN can be used to simulate clusters 
with multiple starbursting events. As a demon- 
stration of this, we present the twin cluster h 
and X Persei (NGC 869, NGC 8 84), simulate d 
based on measurements made by iBraed ( 20041 ). 



Our V-band simulation is shown in Figure 1171 
Here, we applied the following inputs to our sim- 
ulation: \og(ageh) — 7.05, log(age x ) = 7.09, 
M h = 5.5 x 1O 3 M , M x = 4.3 x 10 3 M Q . For both 
clusters we assumed solar metallicity, Ay = 1.6 
mag, Rv = 3.1, and a distance of 2 kpc (dis- 
tance modulus 11.5 mag). The spatial distribu- 
tion of the twin cluster is described by rt = 9.6' 
and r c — 7.01' for h Persei, and rt = 9.6' and 
r c = 8.86' for \ Persei. A Kroupa IMF has been 
used for both clusters, with ot\ = 0.3, «2 = 1-3 
and q?3 = 2.3. A field star population has not 
been applied in this image, making it easier to see 
the extent of the two clusters in the simulation. 

Because we use a King model for the spatial 
extent, it would not be possible to use our code 
to model complex OB associations. However, it 
might be possible to contruct a reasonable OB as- 
sociation using multiple bursts with differing radii 
and age such as demonstrated for h and x Persei. 



(b) 2MASS Image 




Fig. 16. — MASSCLEAN simulated image of the young 
Milky Way cluster, NGC 3603 is given in (a). The 
real image from 2MASS is given in (b). Both images 
represent a 4.4' x 4.4' J Band Image 



Fig. 17. — MASSCLEAN simulated image of the double 
Milky Way cluster, h and X Persei (NGC 869, NGC 
884). The image represents a 30' x 30' field in the V 
band. No field star simulation is included to allow the 
separate cluster radii to be better viewed. 
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3.5. Image simulations of extragalactic 
clusters 



a future version of MASSCLEAN since the individual 
stellar masses are tracked with time. 



Because the linear scale of the cluster (in par- 
sees) and distance to the cluster is determined by 
the user, one can use MASSCLEAN to simulate a clus- 
ter at any distance, even extragalactic stellar clus- 
ters. We have created a series of simulations made 
using the same initial inputs as in Figure fTSt a) and 
Figure [HJa) of Westerlund 1 and NGC 3603, re- 
spectively Only now we have simulated these clus- 
ters over a range of ages \og(age/ years) — [6,8] 
and placed the clusters at the distance of M31. 
The simulation is designed to represent the depth 
and resolution (0.14") expected from the Hubble 
Space Telescope using the Advanced Camera for 
Surveys instrument. Both V-band and I-band sim- 
ulated images are presented. A model for the field 
stars has also been applied. We generated field 
stars close to the ones observed in real cluster im- 
ages in M31 and to keep the cluster distinguishable 
from the stellar field for at least 10 million years. 
Beside each image we provide the current view of 
the CMD for that cluster. All CMDs are scaled to 
the same magnitude limits. 

The images given in Figures US] through [2D 
have been rendered from the original MASSCLEAN 
generated FITS files to pdf. Thus considerable 
dynamic range of the magnitudes originally con- 
tained in those FITS images have been lost. We 
have also been forced to select a single set of mini- 
mum and maximum brightness levels and a single 
slope greyscale solution to view the clusters. This 
allowed us to show the relative change in the clus- 
ter's appearance with age, though the levels may 
not be optimal for all age periods. This has all 
lead to highly degraded images as given in this 
paper. However, one can recover the full dynamic 
range of the magnitudes given by the package by 
viewing the FITS files included in the electronic, 
on-line version of this paper. 

In Figures [18] through [21] we chose a static 
model for the spatial distribution of stars ((keep- 
ing the stars in the same location) in order to 
allow the viewer to easily compare the bright- 
ness of stars at different ages. Naturally, clus- 
ters are prone to many kinematic disruptions, in- 
ternal and external. Kinematics and disruption 
mechanisms, such as the N-b ody code STARLAB 
( Portegies Zwart et al.ll2001 ). could be included in 



4. Applications of MASSCLEAN 

MASSCLEAN is not meant to replace previous, 
widely used SSP models. In fact, in applications 
where spatial resolution is not available and spec- 
tral information is, MASSCLEAN will be an inferior 
choice. However, we have identified a few astro- 
nomical questions of interest to us that drove us 
to create an image simulation routine, based on 
the tenants of the modern SSP models. We will 
discuss just a few of these below. 

4.1. Selection effects in stellar cluster 
searches 

Researchers have been developing systematic 
surveys of large scale Galactic databases in the in- 
frared to identify new, previously unknown stellar 
clusters. Some of the earliest and most success- 
ful of those investigations were c ompleted using 
the 2MASS (ISkrutskie et al.ll2006h surv ey, first by 
iDutra et al1l2003l and iBica et al.ll2003l More re- 
cen tly an extensive study of the 2MASS database 
by iFroebrich et al. 2007 has yielded more than 
1000 new stellar cluster candidates. The Spitzcr 
Spac e Telescope GLIMPS E galactic plane sur- 
vey ([Benjamin et all l2003h and U KIDSS Galac- 
tic plane survey (jLucas et al ] |2007l) have recently 
gone public leadin g to additional cand idate stellar 
cluster detections ( Mercer et alj|2005l) . 

In uncovering the galactic open and globu- 
lar cluster population, researchers have become 
well aware of the problem of false positives. A 
sight-line with slightly reduced average extinction 
will peer deeper into the inner galaxy causing an 
observed stellar density enhancement, mimicking 
a clu ster ( Cotera fc Leistr al 120031 IFroebrich et al 



2007). But false positives are not the only errors 



important in obtaining a complete picture of the 
star clusters of our galaxy. There also lies the error 
of the second kind: the false negative. 

As an example, it has been proposed that very 
massive clusters (mass > 10 5 A/ Q ) do not exist in 
the Milky Way Galaxy, either because the Milky 
Way is not capable of creating very massive clus- 
ters or the power law of the cluster m ass function 
turns over (jGieles et al.l l2006a|2006bl ). Since we 
do not see clusters of this high a mass within our 
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Fig. 18. — Westerlund 1 in M31 (V Band), simulated images and color-magnitude diagrams. 
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Fig. 19. — Westerlund 1 in M31 (I Band), simulated images and color-magnitude diagrams. 
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Fig. 20. — NGC 3603 in M31 (V Band), simulated images and color-magnitude diagrams. 
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Fig. 21. — NGC 3603 in M31 (I Band), simulated images and color-magnitude diagrams. 
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galaxy, do we accept that they do not exist? And 
what if they do exist but are not detectable by 
present surveys? Then we have made the error 
of failing to observe an object when in fact there 
is one. Such an error is called a false negative. 
The extent of false negatives is impossible to de- 
rive directly from the clusters which are detected 
nor can it be estimated analytically. An image 
simulation package such as MASSCLEAN, when used 
with a realistic simulation for the resolution and 
sensitivity of current Milky Way surveys, can be 
used to derive selection effects and expected yields 
from present search methods. 

4.2. Constraining field-star contamination 

Beyond false negatives, efforts to identify and 
characterize even optically observed galactic clus- 
ters c an suffer greatly from field star contamina- 



tion (jBonatto et al.l 120061) . particularly for low 



contrast clusters. Bonatto fc Bica 2007 estimate 



that 10-20 percent of radii of open clusters may be 
underestimated due to the confusion brought on 
by high field star densities. Moreover, with time 
the dynamical evolution of these clusters, leading 
to increased mass segregation, further reduces the 
surface brightness at large radii. Imaging simu- 
lations can begin to address the extent of these 
biases in various fields, in different bands and as a 
function of galactic latitude, stellar field properties 
and cluster mass segregation. 

4.3. MASSCLEAN as part of a Monte Carlo 
analysis 

MASSCLEAN determines integrated luminosities 
and colors, such as those shown in Figure [SJ based 
on the pure summing up of the luminosity of the 
individual objects in the simulation as individual 
objects are being tracked in the simulation. This is 
quite different from the standard SSP codes which 
rely on a probabilistic mass function of stars to 
describe luminosities. In the later method, there 
is no way to estimate the 'stochastic fluctuations' 
expected with a finite number of stars in a real 
cluster and the range of obse rved integrated colors 
which will naturally result (jCervino fc Luridiana 
2006) . The mean values derived from probabilistic 
models are correct on average, but they may not 
represent the actual values observed, particularly 
in clusters which lie b elow the Lowest Luminosity 
Limit as described by ICervino fc Luridianal 20041 . 



MASSCLEAN, when used as part of a Monte Carlo 
analysis, can be used to characterize the expected 
distribution of integrated broad-band colors as a 
function of cluster mass. We demonstrate this al- 
ready in §3.1. 

5. Conclusion 

Our MASSCLEAN package has been introduced. 
We provide a few first order checks that our simu- 
lations are consistent with other, well accepted, 
modern SSP codes. We have demonstrated its 
features simulating various clusters over a mass 
range of a few xlO 3 — 10 5 M Q and showing sim- 
ulated CMDs and single-band images. Although 
the package can simulate low mass clusters, we 
emphasize its use for massive clusters. In order to 
determine the mass, age and distance of a clus- 
ter, observational astronomers typically rely on 
the brightest stars. Perhaps no more than 10% of 
the cluster is visible in typical images, particularly 
when extinction is high. Our approach comple- 
ments the efforts of observers to take full advan- 
tage of all information obtained in their images 
and by simulating all of the stars (not only the 
visible ones). In a forthcoming paper, we will de- 
velop and present methods to derive the goodness 
of fit for a simulated versus a real cluster using 
cumulative distribution functions. 

We are grateful to suggestions made to an 
early draft of this work by Rupali Chandar, Bruce 
Elmegreen, Deidre Hunter and Brad Whitmore. 
Their ideas lead to significant improvements in 
the presentation. This material is based upon 
work supported by the National Science Founda- 
tion under Grant No. 0607497 to the University of 
Cincinnati. BP was supported by a graduate sum- 
mer fellowship from the University of Cincinnati's 
Research Council. 
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